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ABSTRACT 


Current cloud microphysical schemes used in cloud and mesoscale models range from simple 
one-moment to multi-moment, multi-class to explicit bin schemes. This study details the 
benefits of adding a 4* ice class (frozen drops/hail) to an already improved single-moment 3- 
class ice (cloud ice, snow, graupel) bulk microphysics scheme developed for the Goddard 
Cumulus Ensemble model. Besides the addition and modification of several hail processes 
from a bulk 3-class hail scheme, further modifications were made to the 3 -ice processes, 
including allowing greater ice supersaturation and mitigating spurious evaporation/sublimation 
in the saturation adjustment scheme, allowing graupehhail to transition to snow via vapor growth 
and hail to transition to graupel via riming, wet graupel to become hail and the inclusion of a rain 
evaporation correction and vapor diffusivity factor. The improved 3 -ice snow/graupel size- 
mapping schemes were adjusted to be more stable at higher mixing rations and to increase the 
aggregation effect for snow. A snow density mapping was also added. 

The new scheme was applied to an intense continental squall line and a moderate, loosely- 
organized continental case using three different hail intercepts. Peak simulated reflectivities 
agree well with radar for both the intense and moderate case and were superior to earlier 3-ice 
versions when using a moderate and large intercept for hail, respectively. Simulated reflectivity 
distributions versus height were also improved versus radar in both cases compared to earlier 3 - 
ice versions. The bin-based rain evaporation correction affected the squall line more but overall 
agreement in reflectivity distributions was unchanged. The new scheme also improved the 
simulated surface rain rate histograms. 


2 


71 

72 

73 

74 

75 

76 

77 

78 

79 

80 

81 

82 

83 

84 

85 

86 

87 

88 

89 

90 

91 

92 

93 

94 


1. Introduction 


Atmospheric cloud modeling has benefited immensely from the continued improvement in 
computational power. Simulations using explicit spectral bin microphysics (SBM) with large 3D 
domains in mesoscale models like WRF (the Weather Research and Forecasting model, 
Michalakes et al. 2004; Skamarock et al. 2008) can now be performed (Iguchi et al. 2012a, b). 
In addition to higher resolution (e.g., Khairoutdinov and Randall 2006) and the advent of MMFs 
(multi-scale modeling frameworks, Randall et al. 2003, Tao et al. 2009), cloud-resolving 
simulations have also benefited in the form of ever more sophisticated microphysics. Simple 
bulk liquid (e.g., Kessler 1969) and ice schemes (e.g., Wisner et al. 1972) with only a few 
categories have grown into multiple ice categories (e.g., Straka and Mansell 2005), two moments 
(e.g., Ferrier 1994; Reisner et al. 1998; Morrison et al. 2009) and higher (Milbrandt and Yau 
2005b), and highly detailed SBM (Ovtchinnikov and Kogan 2000; Khain et al. 1999; 2000; 
2004). Detailed bin forms originated a while ago but are only now becoming practical, having 
previously been limited to either ID (Young 1974; Scott and Hobbs 1977), 2D (Takahashi 1976; 
Hall 1980; Reisin et al. 1996; Khain and Sednev 1996) or without ice (Kogan 1991). 

Though the ability to use SBM with ever increasing 3D domains is becoming a reality, these 
types of simulations are still not common, and there is a trade-off versus bulk microphysics 
schemes (BMSs), which assume hydrometeor distributions follow a prescribed form (typically 
exponential or gamma). BMSs are much faster and require a lot less memory. For certain 
applications (e.g., very large domains or long-term simulations), computational resources are 
often not sufficient for SBM, which themselves are still not perfect (Li et al. 2010); these 
resources could also be applied to better resolution, which is another important consideration 
when it comes to realistically simulating convective entrainment and overturning (Bryan et al. 
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2003). BMSs are typically invoked using either one-moment (IM, only mass is predicted) or 
two-moments (2M, both mass and number concentration are predicted); these schemes have also 
seen numerous advancements and improvements in recent years. Numerous modeling studies 
and BMSs were made or based on the IM 3-class ice (3ICE) schemes of Lin et al. (1983) and 
Rutledge and Hobbs (1983, 1984) developed in the early 1980’s. These schemes were the 
wor kh orses of cloud microphysics for many years and are still used in some form by many 
schemes today, especially for NWP. However, they do have their biases (Lang et al. 2007, 2011; 
hereafter L2007, L2011) and are susceptible to thresholding phenomena (Rutledge and Hobbs 
1984). They use a priori hydro meteor classes of cloud ice, snow and either graupel or hail and 
transfer hydrometeors from one class to another conditional upon specified thresholds; this can 
result in abrupt and unnatural behavior and diverging solutions depending on if conditions are 
met. An innovative approach was recently developed by Morrison and Grabowski (2008) based 
on the concepts of Heymsfield (1982) and Hashino and Tripoli (2007) whereby the amounts of 
mass acquired by riming and deposition are predicted separately. This allows for the history of 
the riming fraction to be accounted for and results in a spectrum of particle densities with 
smooth, natural transitions from ice to snow and snow to graupel. Lin and Colle (2011) included 
the effects of partially rimed particles using a diagnostic riming intensity as well as functional 
forms of the mass-, area- and velocity-diameter relationships. In general, the representation of 
cloud microphysical processes is constantly improving as more and more schemes are including 
aerosols and the representation of ice processes continues to improve (Muhlbauer et al. 2013). 

Two-moment BMSs have become increasingly popular (Ziegler et al. 1985; Terrier 1994; 
Reisner et al. 1998; Meyers et al. 1997; Milbrandt and Yau 2005a; Morrison et al. 2005, 2009; 
Seifert and Beheng 2006; Thompson et al. 2008; Lim and Hong 2010) and offer a good 
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compromise between the extreme eost of SBM and the restrietion of IM BMSs. They allow an 
extra degree of freedom in defining the hydro meteor distributions compared to IM sehemes and 
ean aceount for size-sorting (Milbrandt and Yau 2005a) and aggregation as well as aerosol 
effects (Lim and Hong 2010). Two-moment schemes are also superior to IM in terms of rain 
evaporation. A fixed rain intereept in IM schemes tends to produee exeessive rain evaporation, 
namely in the stratiform region; this ean be alleviated with 2Ms, which can lead to better 
conveetive and stratiform rainfall struetures (Morrison et al. 2009). 

However, 2M sehemes as well as SBM require observed eloud eondensation nuelei (CCN) 
and/or iee nuelei (IN) profiles to aetivate cloud water and/or ice partieles, and despite their 
potential (Terrier et al. 1995; Milbrandt and Yau 2005a; Morrison et al. 2009), there are enough 
uneertainties and nonlinearities that more advanced/sophistieated sehemes do not always perform 
better than simpler IM bulk schemes (Wang et al. 2009; Varble et al. 2011). Besides the extra 
degrees of freedom that are required to behave in a realistic manner, larger errors in a more 
dominant proeess ean overwhelm potential gains elsewhere. WRF contains a variety of 
mierophysies packages, including IM, 2M and schemes that are a mixture of both, easily 
allowing eomparisons between the various sehemes. Van Weverberg et al. (2013) used the 
Advaneed Research WRF to evaluate MCS simulations using three different BMSs over the 
tropieal western Paeifie and found that although different, the results from the 2M sehemes were 
not superior to those using IM; they found the most crueial element was the fall speeds of frozen 
particles. Using WRF, Han et al. (2013) found that the 2M Thompson seheme (Thompson et al. 
2008) had the best radar reflectivities for a winter storm but the fall veloeities from the Goddard 
scheme (Tao and Simpson 1993; Tao et al. 2003; L2007) agreed best with vertieal profiler 
observations. Powell et al. (2012) used vertically-pointing millimeter radar observations to 
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evaluate WRF simulations of MCS anvils near Niamey, Niger; they found the IM Goddard and 
WDM6 (Lim and Hong 2010, which has IM ice and 2M liquid) schemes actually produced more 
realistic anvils than did the 2M ice schemes. Systematic biases in cloud-resolving model (CRM) 
BMSs were first identified when CRM simulations were used for satellite retrievals. Radiative 
transfer models applied to CRM-simulated cloud fields revealed distributions that contained 
excessive scattering signatures that were not representative of actual observed distributions 
(Panegrossi et al. 1998; Bauer 2001; Olson et al. 2006). In addition to the excessive scattering, 
several studies found excessively high reflectivities in the upper troposphere of CRM simulations 
with excessive amounts and/or sizes of precipitation-sized ice produced by the BMSs as the 
primary reason (L2007, L201 1 ; Blossey et al. 2007; Zhou et al. 2007; Li et al. 2008; Matsui et al. 
2009). Typically, the problem is associated with graupel (L2007; Li et al. 2008; Matsui et al. 
2009). Though somewhat confirmed by Varble et al. (2011), they noted that excessively high 
reflectivities can also be due to snow in a 2M scheme (Morrison et al. 2009) and that reasonable 
reflectivities can be obtained using nonspherical graupel particles with variable density and 
mass-diameter relationships with IM. At any rate, there are enough inherent biases and room for 
improvement that microphysics schemes in general require and continue to undergo refinement. 

This study is a follow on to L2007 and L2011 and details the continued improvement and 
enhancement of the Goddard IM BMS used in the Goddard Cumulus Ensemble model (GCE), a 
version of which (E2007) is one of the microphysics packages available in WRF. This scheme 
has been evaluated in WRF and despite some continued biases (E2007, E201 1) found to compare 
quite well relative to other WRF schemes with regard to snow fallspeeds (Han et al. 2013) and 
MCS anvils. The lineage is built upon the Rutledge and Hobbs (1983, 1984) 3ICE-graupel 
version. Despite the improvements already made to the scheme (E2007, L201 1), it still contains 
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some unrealistic aspects with regard to reflectivity structure and, without hail, lacks the ability to 
simulate more intense radar echoes. Hail allows for the simulation of much more intense radar 
echoes and much higher fall velocities than graupel and can thus to flrst order expand the scheme 
to cover a far wider range of conditions than would be possible with say improving the graupel 
category to 2M. This study details the addition of a 4* ice category, which encompasses the 
spectrum of particles from smaller frozen drops to larger hail stones that have a high density 
(-0.9 g cm‘^) as a result of being or having a coating at or near liquid at some point in their 
history, as well as further reflnements to the scheme that result in an improved 4-class ice (4ICE) 
version of the Goddard IM BMS. Despite the beneflts of a smooth transition in particle 
characteristics (e.g., density) with a single prognostic rimed-ice category, higher and lower- 
density particles cannot coexist in such a scheme without a separate hail category, which may be 
necessary to simulate for example a narrow hail shaft (Milbrandt and Morrison 2013). The new 
Goddard 4ICE scheme is tested for two cases, an intense midlatitude squall line observed during 
MC3E and a more moderate convective case from TRMM LBA (the Tropical Rainfall 
Measuring Mission Earge-Scale Biosphere -Atmosphere Experiment in Amazonia). The model 
results are evaluated using radar reflectivity contoured frequency with altitude diagrams 
(CEADs, Yuter and Houze 1995). Validation via comparison with in situ aircraft data can 
provide a very detailed look at the performance of microphysical schemes (e.g., Molthan and 
Colic 2012); however, such data are limited and difficult to compare against (if even available) 
when it comes to convective cores and are essentially unavailable when it comes to intense 
convective cores. Another approach has been to compare modeled versus observed radiances, 
often radar reflectivity, statistically in the form of CEADs, the primary approach adopted in 
L2007 and L20 1 1 . CEADs, which are essentially PDEs sampled at discrete levels through the 
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depth of a storm stacked in the vertical and then contoured, were first used primarily to 
characterize observations (e.g., Yuter and Houze 1995). Lin (1999) and Lang et al. (2003) 
constructed CFADs of model data, but Smedsmo et al. (2005), Eitzen and Xu (2005) and Braun 
(2006) were the first to use radar CFADs to actually evaluate CRM performance. These were 
quickly followed by several other studies (Blossey et al. 2007; L2007; Rogers et al. 2007; Zhou 
et al. 2007; Li et al. 2008), establishing thi s method (or a close variation thereof) as a standard 
way of evaluating CRM-type simulations (L2011; Varble et al. 2011; Iguchi et al. 2012a; Powell 
et al. 2012; Guy et al. 2013; Han et al. 2013; Van Weverberg et al. 2013; Wu et al. 2013). 
Though there are other ways to evaluate the performance of CRMs and their microphysics 
schemes, having to match the observed radar reflectivity distributions is a more stringent test 
than mean quantities (Powell et al. 2012). Radar observations are more readily available than in 
situ observations, especially when it comes to convection, and the resolution of the data is 
comparable to CRM grids and much better than satellite observations. 

The main objectives of thi s study are to allow the improved IM Goddard BMS (L2007; L201 1) 
to simulate intense convection via the addition of hail and to further improve upon the model’s 
overall performance via the reduction of biases in the synthetic radar structure and reflectivity 
distributions. The paper is organized as follows. Section 2 describes the Goddard CRM, the 
changes to the Goddard microphysics, the two case studies, and the numerical experiments. 
Section 3 presents the simulation results and their validation using radar observations as well as 
surface rain intensities. The summary and conclusions are given in section 4. 

2. Simulation setup and cases 
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a. The Goddard Cumulus Ensemble model (GCE) 


The new Goddard IM 4ICE BMS is evaluated using the 3D GCE. The GCE was deseribed in 
Tao and Simpson (1993) and Tao et al. (2003) and more recently in Tao et al. (2013). The 
model configuration closely follows that used in L2011. The GCE has a one-and-a-half order 
sub-grid scale turbulence scheme (Soong and Ogura 1980), parameterizations for shortwave 
(Chou and Suarez 1999), longwave (Chou and Kouvaris 1991; Chou et al. 1995, 1999; Kratz et 
al. 1998) and cloud optical properties (Sui et al. 1998; Fu and Liou 1993), positive definite 
advection (Smolarkiewicz 1983, 1984; Smolarkiewicz and Grabowski 1990), and options for 
anelastic (Ogura and Phillips 1962) or compressible flow (Klemp and Wilhelmson 1978). The 
GCE has several microphysics options, but the primary BMS is the Rutledge and Hobbs (1983, 
1984)-based 3ICE-graupel scheme (i.e., cloud water, rainwater, cloud ice, snow and graupel), 
which has been improved to reduce unrealistically large amounts of graupel (E2007) and 40 dBZ 
echoes above the Ifeezing level (L2011) and modified to introduce ice nuclei concentrations into 
the Bergeron parameterization (Zeng et al. 2008, 2009). 

b. Addition of hail processes and other microphysics improvements 

Prior improvements reduced excessive amounts of graupel (E2007) and excessive penetrations 
of 40-dBZ echoes above the Ifeezing level (E201 1), alleviating some of the biases in the original 
Goddard IM 3ICE-graupel scheme. The improved versions, however, have two artifacts in their 
simulated reflectivity structure. Time-height cross sections of peak reflectivities show a band of 
elevated reflectivity maxima above the freezing level separated by corresponding local minima 
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just above the freezing level (e.g., Figs. 2c-d, 6d) as opposed to observations, which typically 
show a steady monotonic decrease with height (e.g., Figs. 2a, 6a). The snow/graupel size 
mapping implemented in L2011 tried to compensate for this and still produce reasonable peak 
reflectivities by rapidly increasing the particle sizes (especially snow) at moderate to high mixing 
ratios (Figs. Ic-d). This can lead to spurious artifacts in the reflectivity distributions near the 
melting layer for stronger cases (Fig. 4d). To address these issues and allow the scheme the 
ability to simulate more intense convection, a 4* ice class (frozen drops/hail) was added and the 
scheme further refined to produce an improved Goddard 4ICE scheme. 

Table 1 lists the changes in relation to previous versions of the scheme, which are detailed 
below. Hail processes based on Lin et al. (1983) were added with some modification. Hail 
riming, accretion of rain, deposition/sublimation, melting, shedding and wet growth processes 
were left unchanged. Analogous to L2007 for graupel, the collection of other ice particles under 
dry growth conditions (dry collection) was eliminated for hail to prevent an excessive buildup of 
hail as collection efficiencies should be minimal; but, hail near wet growth conditions is 
expected to be close to water coated and thus efficiently collect other ice particles. Hail within 
95% of wet growth is thus allowed to collect other ice particles. Graupel is medium density 
(-0.3 to 0.5 g cm'^) and mainly a riming product; frozen drops, however, are high-density (-0.9 g 
cm‘^) like that of solid ice or hail. Therefore, the new scheme differentiates hail initiation from 
graupel by treating the product of any process that freezes rain as hail. Five new hail processes 
were also added: wet hail accretion of graupel, hail rime splintering, wet graupel conversion to 
hail, hail conversion to snow due to depositional growth (also added for graupel) to allow 
haiFgraupel particles that undergo significant deposition at cold temperatures in the absence of 
liquid water to transition to snow, and hail conversion to graupel via riming, which tries to 
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account for lower hail bulk densities due to rime accumulation under non wet growth conditions. 
Milbrandt and Morrison (2013) demonstrated how graupel densities can sharply decrease at 
colder temperatures aloft when accounting for variable rime density. A similar analogy is 
applied here with regard to the hail category, which too can be rimed. The latter two processes 
are “apparent” and try to overcome some of the limitations associated with having fixed particle 
categories by providing additional pathways for particles to change Ifom one category to another. 

In addition to the hail processes, further modifications were made to the improved SICE 
processes. Snow autoconversion was strengthened by adjusting the timescale, threshold and 
efficiency in the Kessler-type formulation. The original formulation is quite weak, but aircraft 
observations of small ice particle distributions suggest that although diffusion dominates at 
colder temperatures, there is evidence of aggregation (Field 1999). Though autoconversion was 
strengthened, diffusion still remains the dominant process. The maximum snow collection of 
cloud ice efficiency was increased to reflect the fact that ice particles with diameters greater than 
about 200 microns are efficient collectors. A water vapor diffusivity correction factor (Byers 
1965) was applied to all related processes, effectively replacing a more simplified function of 
pressure and temperature. The formulation for the depositional growth of cloud ice into snow 
(i.e., Psfi) invokes the time step over a time scale, so that tiny amounts of snow can form even 
when cloud ice is quite small (i.e., the time scale to become snow is large). Therefore, an 
arbitrary small threshold was introduced before activating this term^ With ice supersaturations 
commonly observed on the order of tens of percent (Jensen et al. 2001; Stith et al. 2002; Garrett 
et al. 2005), the sequential saturation adjustment scheme (Tao et al. 2003) was further relaxed 
from the 10% constraint in L2011 to a maximum of 20%. The deleterious effects of cloud 

* A value of 1 .e-5 g/g scaled by the surface air density over the level air density was used to inhibit this artificial 
snow production mainly because it has a noticeable effect on longer-term (i.e., MMF) simulations. 
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boundaries advecting through an Eulerian grid (e.g., spurious evaporation) have been previously 
noted [Klaassen and Clark 1985; Grabowski 1989; Grabowski and Smolarkiewicz 1990; Stevens 
et al. 1996; Grabowski and Morrison 2008; Reisner (personal eommunication)]. To reduce these 
effects, the saturation adjustment was further modified to restrict cloud evaporation and 
sublimation to subsidence areas. The time-scale for ice sublimation can be appreciable, allowing 
even smaller ice particles to exist in subsaturated conditions (e.g., Garrett et al. 2005); so, the 
saturation adjustment was relaxed to allow cloud ice to persist under ice subsaturated conditions. 

With the addition of hail, the size-mapping scheme for snow and graupel introduced in L2011 
was adjusted. The rate at which the characteristic size (i.e., inverse of the slope parameter) 
increases with mixing ratio, especially of snow, was lowered (see Fig. 1). As will be shown, 
peak reflectivities are now determined largely by hail and do not require large graupel or snow 
particles to generate these higher values. In addition, the aggregation effect, especially for snow, 
was increased to allow particle sizes to grow more rapidly as temperatures rise from -25 C, and 
an associated snow density mapping (Fig. Ig) was introduced as a function of snow size (i.e., 
Brandes et al. 2007). Graupel density was also divided into low and moderate based on a simple 
mixing ratio threshold. Finally, to address the problem of excessive rain evaporation due to a 
fixed rain intercept, a rain evaporation correction was adopted based on the results from an 
explicit bin microphysics scheme (Fi et al. 2009). This correction was made “physicaf’ by 
lowering the rain intercept (i.e., increasing the mean raindrop diameter) locally for each rain grid 
cell until the new evaporation rate matched the correction factor. The scaling was capped such 
that rain size never increases with decreasing rain mixing ratio. 

For a more detailed description of the new scheme, including the new hail budget equation, 
rain evaporation correction and snow mapping, please see the Appendix. 
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c. Convective case studies 


1) AN INTENSE MIDLATITUDE SQUALL LINE: MC3E 

The Midlatitude Continental Convective Clouds Experiment (MC3E) was conducted in and 
around central Oklahoma (OK) from April to May 2011 as a collaborative effort between the 
DOE (Department of Energy) ARM (Atmospheric Radiation Measurement) and NASA GPM 
(Global Precipitation Measurement) programs. The 20 May 2011 case featured a deep, upper- 
level low over the central Great Basin moving into the central and southern Rockies before 
lifting into the central and northern Plains. At the surface, low pressure in southeastern Colorado 
drew warm, moist air up through the southern Plains to a warm front oriented E-W across 
Kansas, while a dry line extended from the Texas/Oklahoma Panhandle down through the 
Concho Valley. The result was a series of convective lines that formed over the Great Plains and 
propagated eastward toward the Mississippi Valley. The most intense squall line to pass through 
the MC3E sounding network was the result of convection that had developed over south-central 
Kansas (KS) and north-central OK within the network merging with the northern end of a long 
convective line that had formed along the dry line and extended from southwestern OK down to 
the Big Bend. The northern portion of this longer line entered the MC3E sounding network 
around 07 UTC 20 May and by 09 UTC had consolidated with the convection near the KS-OK 
border to form a more intense convective segment with a well-defined trailing stratiform region 
that then propagated through the network between 09 and 12 UTC. By 13 UTC, the convective 
leading edge had exited the network, leaving the network dominated by a large area of stratiform 
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rain. Peak reflectivities within the network exceeded 60 dBZ over a depth of several km from 2 
km above and to 3 km below the freezing level, a strong sign of signiflcant hail (May and 
Keenan 2005; Lerach et al. 2010), while 40 and 50 dBZ echoes reached upwards of 16 and 10 
km, respectively, placing this case within the top 0.005% in intensity (i.e., 1 feature in 20,000 is 
as strong) with regard to 40 dBZ echo penetrations based on a TRMM observed precipitation 
feature climatology (Zipser et al. 2006). With 45 dBZ echoes reaching more than 10 km above 
the freezing level, a further indication of the presence of hail (Waldvogel et al. 1979); this case is 
well suited to test and evaluate the new Goddard IM 4ICE scheme with hail. 

As in previous GCE modeling studies (e.g., Zeng et al. 2007, 2008; E2011), the 3D GCE was 
driven by large-scale forcing data obtained from a variational analysis approach (Zhang et al. 
2001), in this case, from the MC3E sounding network. Model experiments were run for 4 days 
starting at 0000 UTC 17 May and ending at 0000 UTC 21 May using 1-km horizontal grid 
spacing and a 256 x 256 km horizontal domain (similar in size to the sounding network) with a 
stretched vertical grid having 76 levels and a top near 27 km. A 2-second time step was used, 
and surface sensible and latent heat fluxes were imposed from the variational analysis. 

2) A MODERATE TROPICAE CONTINENTAL CASE; TRMM LBA 

The 23 Eebruary 1999 case was previously presented in L2007 and L2011. It was 
characteristic of the widespread, weaker monsoon-like convection observed within the westerly 
wind regime (Cifelli et al. 2002; Rickenbach et al. 2002) during TRMM LBA. It falls within 
about the top 1% of TRMM precipitation features (i.e., 1 feature in 100 is as strong) in terms of 
40 dBZ echo height intensity (Zipser et al. 2006) and provides a good contrast to the more 
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intense 20 May 2011 MC3E case. On thi s day, daytime heating triggered widespread convection 
that loosely organized into southeast-northwest bands. A long thin convective band developed 
by 2000 UTC and by 2100 UTC was already decaying. Dual-Doppler observations were 
collected for the northern portion of this line where 40 dBZ echoes reached to about 7 km. 
Please see L2007 and L201 1 for more details. 

The current model setup closely follows L2011 and uses the same horizontally homogenous 
initial conditions based on the 1200 UTC morning sounding taken at Rebio Jam having weak 
low-level northwesterly flow and a 500 m mixed-layer CAPE of 934 J kg’', cyclic lateral 
boundary conditions and convection initiated by imposing time-varying (diurnal) surface fluxes 
based on surface observations collected at two different sites (ABRACOS Hill and Ji Parana). 
The horizontal domain was kept at 128 x 128 km as in L2011, but the horizontal grid resolution 
was improved to 200 m in both directions and the time step reduced to 2 seconds. The stretched 
vertical grid was kept at 70 levels with a top near 23 km. The north-south oriented rectangular 
patch of higher sensible/lower latent heat fluxes (Ji Parana) imposed in L201 1 was lengthened to 
18 X 80 km (solid rectangle in Pig. 8b) and replaced by higher latent/lower sensible heat fluxes 
(ABRACOS) if accumulated rainfall exceeded 3 mm over the patch to allow for some cloud 
feedback. Simulations were mn for 6 hours. 

d. Numerical experiments 

For each case, seven numercial experiments are conducted. Three experiments are made using 
previous verions of the 31CE-graupel scheme: the original (3ice0), L2007 (3icel), and L2011 
(3ice3). Four variations of the new 41CE scheme are tested: with smaller-, medium- and larger- 
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sized hail with the bin rain evaporation correction (4iceb sml, 4iceb med, and 4iceb Irg, 
respectively) and smaller-sized hail without the evaporation correction (4ice sml). Smaller-, 
medium-, and large-sized hail use fixed hail distribution intercepts of 0.0200, 0.0020, and 0.0002 
cm"^, respectively. See Table 2 for a list of the numerical experiments performed for these cases. 

3. Simulation results and validation 

a. MC3E 

Figure 2b shows the 4ICE control simulation (4iceb sml) compared to the observed convective 
line (Fig. 2a) as it passed through north-central OK. Despite the restrictive double cyclic 
boundary conditions, the model captures the general organization and intensity of the system 
with an eastward propagating north-south oriented intense convective leading edge with a 
trailing stratiform region off to the northwest. A vertical east-west cross section through the 
center of the model domain (Fig. 2c) does show an erect intense uni-cellular convective 
structure, but the simulated trailing stratiform region appears to be noticeably narrower. 

Time-height cross sections (Fig. 3) of NEXRAD and simulated peak radar reflectivity^ within 
the MC3E sounding array and model domain are shown from 00 EITC 20 to 00 EITC 21 May 
2011; the observed period 06 to 12 UTC (Eig. 3a) covers the formation of the main convective 
line within the sounding array until the leading edge propagated out of the array. Peak 
reflectivities within this line exceeded 60 dBZ with 50 dBZ echoes reaching 9 km and 40 dBZ 
echoes 15 km. At 06 UTC, the simulations are still quite weak. The model imposed large-scale 

^ Simulated radar reflectivities were calculated from model rain, snow, graupel and hail contents assuming inverse 
exponential size distributions and accounting for all size mappings using the formulation of Smith et al. (1975) and 
Smith (1984). 
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forcing is horizontally uniform and first results in a patchwork of smaller convective cells over 
the domain that require 2-3 hours to respond to the organizing shear and form into a squall line 
(~09 UTC). The SICE graupel runs (Figs. 3b-d) significantly underestimate the peak intensity of 
the observed squall line above the freezing level. The Sicel and 3ice3 runs produce a band of 
elevated reflectivity maxima 3 km above the freezing level whereas the observed reflectivities 
monotonically decrease with height above the freezing level. Graupel with its relatively slow 
fallspeeds (Fig. 4b) is carried aloft by the strong updrafts in the convective cores (Fig. 4a) where 
its mass is maximized well above the freezing level. In contrast, the 4ICE simulations (Figs. 3e- 
h) produce higher reflectivity values just above the freezing level due to higher fall speeds (Fig. 
4e) that keep the peak hail mass nearer the freezing level as well as peak values that decrease 
monotonically with height: both in good agreement with the observations. Profiles of peak 
reflectivity within the model domain over the period 09 to 15 UTC (sampled to match the 
observed squall line structure from 06 to 12 UTC, Fig. 5) show all three 3ICE simulations have a 
pronounced low bias that ranges from about 5 dBZ below the freezing level to as much as 15 
dBZ above the freezing level. The 4ICE simulations show a marked improvement in the bias at 
almost all levels except for 4iceb Irg, which produces excessively large reflectivities (~I5 dBZ) 
near the melting level. The medium hail profile has the smallest overall bias and agrees best 
with the observed. Though not quite as good, the smaller hail runs are significantly improved 
over the 3ICE with a consistent low bias of just 5 dBZ at all levels. 

Figure 6 shows statistical CFADs constructed over these same observed (06 to 12 UTC) and 
model (09 to 15 UTC) time periods and domains. The observed CFAD shows higher 
concentrations broadly ranging from 0 to 40 dBZ below the melting level; aloft a coherent core 
of higher probabilities increases from 10 dBZ near 200 mb to 25 dBZ just above the freezing 
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level. Concentrations of infrequent but more intense echoes extend out to near 65 dBZ at and 
below the freezing level, 50 dBZ at 12 km, and 40 dBZ at 16 km. None of the simulated SICE 
CFADs (Figs. 6b-d) produce reflectivities over 50 dBZ above the freezing level and thus miss 
the stronger echoes in the tail of the observed distribution consistent with Fig. 5. Furthermore, 
the core of highest probabilities is either too broad with too many 40 dBZ echoes due to graupel 
(SiceO), shifted too high with too many 30 dBZ echoes due to snow (Sicel), or shifted too low 
with too many weak echoes (SiceS) compared to the NEXRAD data. In contrast, all of the 4ICE 
simulations produce much better concentrations of both the moderate core echoes and the 
infrequent but intense echoes that arise from hail relative to the observations (Figs. 6e-h). The 
4ICE core distributions are still too broad as a result of too many weak echoes, but their overall 
slope along on the right edge is fairly well aligned with the observed core probabilities. 

The improvements in the 4ICE distributions are reflected in the normalized overlap score 
between the observed and simulated PDFs at each level in the CFADs where unity represents 
perfect overlap and zero no overlap. Figure 7 shows the 4ICE PDFs are consistently better than 
the SICE above the freezing level from 6 to 10 km. Below the melting layer, all of the 
simulations are similar and better than in the mixed and ice phase regions. The agreement 
between the simulations and observations drops off sharply near storm top where simulated radar 
echoes become too weak. This discrepancy was noted in F201 1 and could be due to entrainment 
effects wherein dry air disproportionately sublimates small particles while preserving relatively 
large particles. Cloudsat CFADs show a distinct difference between convective and anvil 
regions (Fuo et al. 2009) with the highest probabilities in the anvil concentrated at the lowest 
reflectivities, indicating mostly smaller particles, before shifting to higher values lower in the 
cloud consistent with accretion and aggregation. In contrast, convective cloud top PDFs are 
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broader, indicating a mix of small and large particles, with a much great proportion of larger 
particles. The current model CFADs resemble the Cloudsat anvil distributions at upper levels. 


b. LBA 


Hail is associated with intense convection and as shown needed to produce stronger radar 
echoes. However, a key objective is for the 4ICE scheme to respond appropriately to the 
environment without the need for excessive tuning. The moderate intensity 23 February case is 
suitable for evaluating the new scheme for a weaker convective environment. A radar CAPPI for 
this case (Fig. 8a) shows the northern end of the transient convective line as it was starting to 
decay. Individual convective cells are loosely aligned with a small stratiform area extending 
northwestwards. The simulated convective leading edge using the 4ICE scheme (4iceb sml, Eig. 
8b) is also cellular in nature and loosely organized into a north-south line with a small stratiform 
area extending northwestward, consistent with the weak southeast-northwest oriented shear on 
this day. A vertical east -west cross section (Eig. 8c) through the center of the domain shows the 
leading edge is somewhat multi-cellular (also see Fig. 10a) with a small, undeveloped stratiform 
area. Following E201 1, a 64 x 64 km subdomain over the northern portion of the simulated line 
(dashed box in Fig. 8b) was used in the analyses with model data averaged to the 1- km resolution 
of the radar analyses. The mean convective fraction within this subdomain over the simulation 
period from 300 to 360 minutes using the 4ICE scheme ranged from 0.46 to 0.49 in close 
proximity to the radar value of around 0.4^. 


^ As in L2007 and L2011, convective fractions were computed based on Rickenbach and Rutledge (1998), a texture 
algorithm applied to radar reflectivity data that largely follows Steiner et al. (1995) to match the radar observations. 
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Figure 9 shows time-height eross seetions of peak refleetivity observed by the NCAR S-pol 
radar within the dual-Doppler analysis domain and simulated within the model subdomains. As 
previously noted (L2007; L2011), the original scheme (3ice0) with its excessive graupel 
produces 40-dBZ echoes that penetrate much higher (over 12 km, Fig. 9b) than was observed (~7 
km, Fig. 9a). Though better, 3icel (Fig. 9c), which eliminates the dry collection of ice/snow by 
graupel (L2007), still produces excessive 40-dBZ echoes penetrations. The 40-dBZ echoes in 
3ice3 (L2011) are greatly reduced aloft and closer to the observed, but 3ice3 results in an 
elevated reflectivity maximum above the freezing level (Fig. 3d) that was not observed. None of 
the 3ICE graupel simulations can reproduce the 45-dBZ echoes immediately above the freezing 
level (~4.9 km) as was observed (Fig. 9a). In contrast, despite the addition of the higher density 
frozen drops-hail ice class, none of the 4ICE simulations (Figs. 9e-h) produces the excessive 40- 
dBZ echo penetrations in 3ice0 and 3icel and all can replicate the observed monotonically 
decreasing reflectivity structure above the melting level though clearly the medium (4iceb med) 
and larger (4iceb Irg) hail peak reflectivities that are too strong around the freezing level. These 
results suggest that hail or frozen drops with their higher fall speeds (Fig. 10c) relative to graupel 
(Fig. 10b) are crucial not detrimental for reproducing the observed core reflectivity structure of 
even moderate convection. Peak reflectivity profiles from the 64 x 64 km model subdomains 
over the final 60 minutes (Fig. 11) confirms the strong (nearly 15 dBZ) to moderate (~5 dBZ) 
over bias in the 3ice0 and 3icel simulations aloft, respectively, as well as the under bias (~8 
dBZ) in 3ice3 near the melting level, which contributes to the elevated maxima. The 4ICE 
simulations with smaller hail clearly perform the best and show almost no bias (less than -4 
dBZ) through nearly the entire depth of the storm. Remarkably, none of the 4ICE runs produce 
the over bias evident in runs 3ice0 and 3icel in the upper part of the storm, and all produce 
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monotonically decreasing profiles with height in agreement with the observations. However, 
quite obviously the medium to larger hail sizes in runs 4iceb med and 4iceb hg are much too 
large, producing over biases of up to ~10 to 15 dBZ around the melting level. These results 
suggest the new 4ICE scheme is quite capable of responding appropriately to the intensity of the 
convective environment and can outperform the 3ICE-graupel scheme in terms of peak 
reflectivities even in a moderate intensity environment. This actually is consistent with 
polarimetric radar and wind profiler evidence that Ifozen drops or hail quite often are present in 
tropical convection (Jameson et al. 1996; May et al. 2001); in situ aircraft data also confirmed 
the presence of Ifozen drops in this case (Stith et al. 2002) 

Simulated and observed radar CFADs (Fig. 12) also show the 4ICE scheme equals or 
outperforms the 3ICE-graupel scheme in terms of overall reflectivity distributions for the 
moderate 23 February case. The highest observed concentrations (Fig. 12a) gradually decrease 
Ifom between ~5 to 20 dBZ at and below the melting level (-d.O km) to — 5 to 15 dBZ near 
storm top with less apparent aggregation than MC3E (Fig. 6a). Fow probabilities of moderately 
strong echoes reach ~50 dBZ at and below the Ifeezing level, 40 dBZ at --7 km, and 30 dBZ at 
12 km. As in MC3E only more pronounced, the original scheme (3ice0, Fig. 12b) produces 
excessive concentrations of 20 to 40 dBZ echoes between the melting level and 10 km, and 
again, though the number of excessive 30 to 40 dBZ echoes is reduced in 3icel (Fig. 12c), the 
core of peak probabilities is ~10 dBZ too high. Core probabilities for 3ice3 are noticeably better 
(Fig. 12d), but the amount of echoes too weak is higher as is the penetration of 20 to 30 dBZ 
echoes near storm top. There is also an unphysical notch in the higher echo distribution just 
above the Ifeezing level. Though hail amounts are small in the 4ICE simulations, they dominate 
the higher reflectivity values; the peak concentrations are nearly invariant (Figs. 12e-h) but are as 
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good as or better than the SICE at every level, aligning well with the observed coneentrations 
along 20 dBZ from the melting level to 12 km. In terms of echoes stronger than 20 dBZ in the 
distribution tail, the smaller hail results (Figs. 12e-f) match the observed frequencies extremely 
well while the medium (Fig. 12g) and larger hail (Fig. 12h) runs are overly intense especially 
near the melting level. And as with MC3E, simulated probabilities collapse below 0 dBZ at 
storm top and are ~10 dBZ weaker than observed. 

The overall level-by-level performance is confirmed by the profiles of normalized PDF 
matching scores (Fig. 13). From the freezing level to 10 km, the 3ice3 and 4ICE PDFs are all in 
excellent agreement with the observed with matching scores on the order of 0.8, much better 
than 3ice0 and 3icel. Above 10 km, the performance of 3ice3 drops off quickly while the 4ICE 
simulations continue to perform well up to 12 km before they begin to deviate from the observed 
echo distributions near storm top. Below the melting level, the 3ICE and 4ICE schemes all 
perform about the same and reasonably well with scores of around 0.75. 

c. Rainfall comparisons and validation 

In addition to radar reflectivity, the 4ICE scheme is validated with regard to surface rain rates. 
Figures 14a,b show the instantaneous surface rain rates associated with the 4iceb sml simulations 
shown in Figs. 2b and 8b for MC3E and FBA, respectively. Surface rain features in the 
convective region are larger, better organized and more intense in MC3E; the stratiform region is 
also larger and more coherent. Simulated surface rain rate histograms (Figs. 14c and d) can be 
compared to observed rain rate histograms derived from the national Doppler radar network for 
MC3E (Fig. 14c) and ground-based radar deployed for the FBA field campaign (Fig. 14d). The 


22 


528 

529 

530 

531 

532 

533 

534 

535 

536 

537 

538 

539 

540 

541 

542 

543 

544 

545 

546 

547 

548 

549 

550 


SICE and 4ICE histograms tend to fall into two distinct clusters, which is more apparent in the 
MC3E case. For MC3E, the SICE simulations significantly underestimate the occurrence of 
more intense rain rates in the tail of the distribution (Fig. 14c); the 4ICE histograms also 
underestimate probabilities of extreme rain rates but are distinctly better than the SICE. For 
LBA, the results are noisier due to a smaller sample size, but the overall results similar: both sets 
of simulations underestimate the proportion of strongest rain rates but with the 4ICE simulations 
closer to the observed histogram than the SICE. These results suggest the 4ICE simulations are 
producing more realistic surface rain rate distributions than the SICE in each environment. 

4. Summary and conclusions 

The improved Goddard SICE IM BMS based on Rutledge and Hobbs (L2011, cloud ice, snow 
and graupel) was modified and hail processes added to produce an improved IM 4ICE BMS 
(cloud ice, snow, graupel, and frozen drops-hail) capable of more realistically simulating the 
radar reflectivity patterns of intense and moderate convection better than previous SICE versions. 
Hail processes taken from the SICE-hail scheme based on Lin et al. (1983) include hail riming, 
accretion of rain, deposition/sublimation, melting, shedding and wet growth. Hail collection of 
other ice species under dry growth conditions was eliminated to prevent the same excessive 
buildup as had occurred with graupel (L2007); however, hail near wet growth is permitted to 
collect other ice particles. Processes that freeze rain now initiate hail not graupel. Five new hail 
processes were added: wet hail accretion of graupel, hail rime splintering, hail (also graupel) 
conversion to snow via depositional growth at colder temperatures, hail conversion to graupel via 
riming under non wet growth conditions, and wet graupel conversion to hail. 
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Besides adding the frozen drops-hail category, snow autoconversion was strengthened based 
on the evidence of aggregation at colder temperatures (Field 1999), the collection efficiency of 
cloud ice by snow increased, a water vapor diffusivity correction factor added, and a small 
threshold introduced to prevent ice deposition growth into snow when ice is small. Maximum 
ice supersaturation was increased to 20%, and cloud evaporation restricted to areas of subsidence 
to mitigate spurious evaporation effects at cloud boundaries. Cloud ice can now persist in ice 
sub-saturated conditions as commonly observed. Snow and graupel size mappings from L201 1 
were adjusted, including an increased aggregation effect for snow. A corresponding snow 
density mapping was added (Brandes et al. 2007), and graupel divided into low and moderate 
densities. Lastly, an SBM-based rain evaporation correction factor (Li et al. 2009) was added. 

The new Goddard 4ICE scheme was used to simulate an intense continental squall line 
observed during MC3E to evaluate its ability to simulate intense convection with significant hail 
as well as a loosely-organized transient line of moderate convection from TRMM LBA to ensure 
the scheme does not over predict less intense convection. For the intense squall line, the 4ICE 
scheme with smaller and medium-sized hail, outperformed prior versions of the 3ICE-graupel 
scheme in terms of peak reflectivities; larger hail produced excessively high values. Without 
hail, the 3ICE-graupel versions could not generate reflectivities over 50 dBZ below freezing, 
vastly underestimating the peak observed reflectivities throughout the depth of the storm with a 
low bias near 15 dBZ at the freezing level, 10 dBZ at midlevels, and 5 to 15 dBZ at upper levels. 
In contrast, the bias was significantly reduced in the 4ICE runs above the freezing level, except 
for the strong over bias of 15 dBZ near the melting level for larger hail. For medium hail the 
bias is less than 5 dBZ at almost all levels except for a positive bias of ~5 dBZ at the freezing 
level. The smaller hail simulations produced a consistent low bias of just 5 dBZ at nearly every 
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level, still a noticeable improvement over the SICE graupel simulations. The 4ICE simulations 
produced radar reflectivity CEADs in better agreement with observations (as reflected in their 
normalized PDE overlapping scores) from 5 to 10 km with more realistic extremes in the 
distribution tails and more realistic reflectivity structures above the freezing level with peak 
reflectivities monotonically decreasing with height as observed versus the SICE graupel 
simulations, which often produced elevated reflectivity maxima. Below the melting level, the 
4ICE runs with smaller hail had peak reflectivities similar to the SICE. Those for medium and 
larger hail were greater due to contributions from melting hail and were closer to the 
observations for medium hail but not for larger where values were excessive. The agreement 
between simulated and observed CEADs below the melting level was similar for all runs. 
Overall, the 4ICE simulation with medium hail performed the best for the intense MCSE case. 

For TRMM LBA, adding frozen drops/hail per se did not necessarily cause unrealistically 
large reflectivity values'^. While medium and larger hail did result in excessive peak 
reflectivities by as much as 10 to 15 dBZ near the melting level, smaller hail had very small 
biases (less than 5 dBZ) that were on average as good or better than the SICE versions at every 
level. Small hail eliminated the over bias at middle and upper levels in the original and L2007 
SICE runs and outperformed the L20I I run above the freezing level by reproducing the observed 
monotonic decrease with height and eliminating the unrealistic elevated reflectivity maxima. 
The 4ICE simulations also produced radar CEADs whose normalized PDF scores were equal to 
or superior to the SICE at all levels. The 4ICE simulations with smaller hail performed the best 
overall for the moderate intensity LBA case. 


"" Trace amounts of mass first appear in the 4* ice class as frozen rain at ~ 10 m/s for the LBA cases and ~5 m/s for 
MC3E, but the sizes are very small and combined with the fact that reflectivity values take a while to grow strong 
shows that the 4ICE scheme is not predisposed to generating larger hail particles nor strong dBZs for every case. 
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Surface rainfall histograms were also used to evaluate the schemes for both cases with similar 
results. In each case the 3ICE simulations significantly underestimated the occurrence of higher 
rain rates compared to observed histograms while the 4ICE histograms had a noticeably higher 
occurrence of stronger rain rates which were closer to but still less than the observed. 

Though the bin rain evaporation correction (Li et al. 2009) did alter the temporal variation of 
peak reflectivities, it had very little effect on either the peak reflectivity profiles or the model 
CFADs in either case despite evidence that it reduced the intensity of the cold pool distribution 
(see Fig. 15). The double cyclic lateral boundaries likely dampened its impact as initial 
excessive evaporation without the correction could over moisten the sub-cloud layer and inhibit 
successive over evaporation thus masking the effects of the correction. The LBA case was in a 
moist environment and of short duration, which could also reduce the impact of the correction. 

The 4ICE scheme with a frozen-drops hail category can simulate more intense radar echoes, 
though choosing a priori the hail intercept for intense or moderate convection is not optimal. A 
size-mapping scheme may alleviate the issue, but ultimately a multi-moment scheme is likely 
needed. 4ICE also replicates the observed monotonic decrease in peak reflectivities with height 
as a result of increasing the range of particle fall speeds to include the higher values associated 
with frozen drops/hail, allowing a greater portion of their mass to remain near the freezing level. 
The scheme also adds pathways (apparent processes) by which particles can move to other 
categories due to their growth mechanisms. Hail can be rimed. If cloud freezes quickly, it 
creates air pockets, which should lower hail density. Just as graupel is assumed to increase in 
density and become hail when reaching wet growth, so too should hail density go down when 
riming rates fall below wet growth. This is the basis for the new Primh hail riming to graupel 
process. For the MC3E control case, Primh increases (decreases) the peak average graupel (hail) 
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content by 25% (5%), decreases the peak echoes by 3 dBZ at and above 12 km, and improves the 
CFAD score by 2 to 4% from 4 to 13 km. The other new set of pathways relates to trace 
amounts of graupel and hail that often persist aloft in the model long after convection ends. 
Initially, sedimentation lowers their mass, but as the mass get small, with a fixed intercept, so 
too do their fallspeeds. The result is tiny amounts of graupel and hail (mean volume diameters of 
a few hundred microns or less) suspended over a broad area where they continue to grow from 
deposition. The new Pvapg and Pvaph vapor conversion to snow processes reduce the area and 
trace amounts of suspended graupel and hail but have minimal impact on the hail cores. For the 
MC3E control case, peak reflectivity profiles are -unchanged by Pvapg/Pvaph below 12.5 km 
(the same is true for Primh below 9 km) as the collection of large amounts of super cooled water 
dominates growth. Peak average snow, graupel and hail amounts inerease -5%, decrease 5%, 
and increase 8%, respectively, when Pvapg/Pvaph are activated; Pvapg/Pvaph improve the 
CFAD score by 2 to 5% from 4 to 1 1 km. The new apparent processes have a small effect on the 
hail cores but allow the model to better address variations in particle density while slightly 
improving the overall echo pattern. 

Allowing a smooth transition in ice particle density (Morrison and Grabowski 2008) and 
fallspeed coefficients (Lin and Code 2011) is a recent and realistic addition to BMSs, though 
comparisons with radar observations for intense cases are needed to evaluate their performance 
and in contrast to having a separate hail category. 

Ultimately, cases should be tested at higher resolution (100 m) to ensure the dynamics are well 
resolved and not contributing to any biases. Further study is also needed to address the 
abundance of weaker model dBZs near storm top to determine if this is an artifact of the radar 
observations or the microphysics. Testing the 4ICE scheme in other environments using remote 
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sensing data for validation (Matsui et al. 2009) is important for systematically identifying and 
eliminating any remaining biases. The ability to match the distribution and peak values of radar 
reflectivities at all levels of a convective system is a fairly stringent test, but radar intensities are 
not a unique solution and can arise from a variety of particle combinations. This is where in situ 
aircraft observations could be very valuable in helping to further constrain the particle 
characteristics. The next iteration of the Goddard scheme is planned to be 2M to ingest and 
include the effects of aerosols and improve limitations inherent in IM. 
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APPENDIX 


a. The hail budget equation 


The hail budget equation for the new 4ICE seheme is given as: 



sacr 


+(l - ^ 2 ) Pi" + Dgacr + Pr acg + Dhacw+ Dhacr + Whaci + Whacs + Whacg + Pglh 


+Phdep - Phsub - Pvaph - Pr imh - Phmlt 


(A1) 


where the first three terms on the RHS of Eq. A1 are the horizontal adveetion, vertieal adveetion 
and diffusion of hail, respectively, while Phfr is the freezing of rain to hail, Piacr cloud ice 
accretion of rain, Prad rain accretion of cloud ice, Psacr snow accretion of rain, Vxacs rain 
accretion of snow, Dgacr graupel accretion of rain, Pr acg rain accretion of graupel, Dhacw 
hail riming, Dhacr hail accretion of rain, Whaci , Whacs , and Whacg wet hail accretion of 
cloud ice, snow and graupel, respectively, Pglh the conversion of graupel to hail, Phdep hail 
deposition, Phsub hail sublimation, Pvaph the conversion of hail to snow via deposition, 
Fvimh hail riming to graupel and Phmlt hail melting. Phfr, Dhacw, Dhacr, Whaci, Whacs, 
Phdep , Phsub, and Phmlt follow the formulations of Ein et al. (1983), while Piacr , Pr act, 
Psacr , Vracs , Dgacr, and Pr acg follow the formulations of Rutledge and Hobbs (1984). 
Whacg follows the Ein et al. (1983) formulation for Whacs but using graupel instead of snow 
parameters. Graupel is assumed to increase in density and become hail upon reaching wet 
growth such that 
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q 

Pg2h = — at the moment when Dgacw + Dgacr > Pgwet 
dt 

where the wet growth of graupel Pgwet is eomputed using the formula for hail wet growth 
from Lin et al. (1983) but with graupel parameters^. 

Just as lower density particles can transition to a higher density class of particles, in the new 
4ICE scheme, the reverse can occur. As such, when hail particles experience riming or 
deposition at colder temperatures, they are transitioned towards graupel and snow, respectively. 
For both processes, an increasing proportion of the mass of hail acquired via riming and 
deposition along with an equal portion of the previous hail mass is transferred to graupel and 
snow, respectively, using a sliding temperature scale with the proportion increasing with 
decreasing temperature. The conversion of hail to graupel via riming is thus formulated as: 

Pr imh = x Dhacw ( A6 ) 

where is given by = 2.0 * (j^f (A7) 

where Tairc is the air temperature in degrees Celcius, tOO is 238.16 degrees K, and tO is 273.16 
degrees K. Similarly the conversion of hail to snow via deposition is formulated as: 

Pvaph = X Phdep when the cloud water mixing ratio qc<i .e-5 g/g (A8) 

where F^^^ uses the same temperature scaling as . The same form of the relation is used for 
the conversion of graupel to snow via deposition Pvapg . Although Hallet-Mossop rime 
splintering (Hallet and Mossop 1974) is not directly part of the hail equation, it does affect the 
hail riming term Dhacw and is computed as: 

^ In nature, larger graupel particles would reach wet growth first. Converting only the mass above the mean volume 
diameter (0.5004) for the MC3E control case reduced the maximum average hail content by over 15%. However, 
peak reflectivities below 9 km were nearly unchanged due likely to high riming rates in the updraft cores but were 
decreased by 3 dBZ or more above 12 km. 
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Pihmh = * Dhacw * 1 000. * Xnsp In ? * Xmsp In t 


where 


?>,„=0.5 

when C <Tairc <-2°C 

II 

o 

when -6° C < Tairc < -4° C 


where the peak number of ice splinters generated per milligram of rime Xnsp\n.t = 370 and the 


mass of each splinter Xmsphit = 4Ae ** grams. Pihmh is first subtracted directly from Dhacw . 
The water vapor diffusivity (Dv) in air was assumed to be a constant (2.26e'^ m^s’') in Rutledge 
and Hobbs (1984); it is parameterized as a function of temperature and pressure in this study 
following Massman (1998): 


T>, - A 


/ \ 1.81 , . 

I P ] f P \ 


v^oy 


5 

\Pj 


where Do is the water vapor diffusivity (2.18e‘^ m^s'^) at To=273.15 K and Po=1013.25 hPa= The 
importance of water vapor diffusivity on the diffusional growth rate of ice crystals as a function 
of pressure and temperature is illustrated in Fig. 9.4 of Rogers and Yau (1989). An adjustment 
factor (Fdwv) is used to adjust the diffusional growth rate of ice crystals based on a constant 
water vapor diffusivity. 


b. Rain evaporation correction 


The rain evaporation correction uses the empirical formula of Li et al. (2009): 
r(^,) = 0.11g“'V0.98 (A9) 
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where r is the ratio of the rain evaporation rate between bulk and bin mierophysies and q,- (g/kg) 
is the rain mixing ratio; r is based on cloud-resolving model simulations using both bulk and bin 
mierophysies and can be used to scale down the bulk rain evaporation rate. In the new 4ICE 
scheme, the correction factor is made “physicaf’ by scaling the rain intercept (i.e., increasing the 
grid local raindrop size) until the bulk rain evaporation rate matches the bin such that: 


where Ftnw is the scaling factor for the rain intercept parameter, 
c. Snow mapping 

The snow mapping scheme maps the snow intercept parameter as a combination of variations in 
temperature and mixing ratio. Variations in mixing ratio are set for two distinct conditions: (1) 
at cold temperatures where aggregation effects are small, sizes are small and only slowly 
increase with increasing mixing ratio and (2) near the melting layer where aggregations effects 
are large, sizes are larger and size increases significantly with increasing mixing ratio. Another 
set of parameters controls how quickly the cold setting variations transform to the warm setting 
variations through the aggregation zone. For both cold and warm regions, an exponent is used to 
control the snow intercept; as this exponent approaches zero, snow sizes relax to those for a fixed 
snow intercept (i.e., larger sizes), and when the exponent approaches one, snow sizes collapse to 
that of a small base size. The formula for the snow size exponent is given by: 


which is then used create the mixing ratio component of the snow intercept scaling factor: 


Ftnw = - 
\r) 



(A10) 



(All) 
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Ftnsq = 


V ^hase J 


(A12) 


where qs\ is the snow water content in g/m^, and Sno\ is a snow water content threshold in 
g/m^ above which snow sizes begin to increase. Snow sizes then continue to increase at an ever 
increasing rate over the next dsn(A g/m^ until reaching the limit is arbitrarily given 

a number close to but less than one. This allows snow sizes to vary ever so slightly (i.e., not be a 
constant size) between snow contents of 0.0 and Sno\ . The parameter settings for cold 
conditions transform through the snow aggregation zone (~ -20 to O'’ C) to those near O'* C as: 


P = Pwarm -(^Pwarm - Pcold ) 


^ Tairc 


Tcold 


(A13) 


where Pwarm is the parameter value near the melting level, Pcold the parameter value for cold 
conditions, and Tcold the air temperature in degrees C for the cold parameter settings. An air 
temperature component for the snow intercept scaling factor, given by 

\ /^exp 


FtnsT = {QXY>{-l.^Tslopes^Tairc)^ 


(A14) 


where Tslopes is the rate of snow intercept change with temperature on a natural logarithm scale 
and Tairc is capped by Tcold , is then combined with the mixing ratio component of the snow 
intercept scaling factor to obtain the total snow intercept scaling factor: 

Ftns = Ftnsq x FtnsT ( A 1 5 ) 

with the condition that the snow size cannot go below a minimum value Dsnowmin- The snow 
intercept mapping is combined with the Brandes et al. (2007) relation between snow density and 

/ ^ S, 0.2995 

median snow volume to get: = 0.001 996 (A16) 
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766 where is snow bulk density and tns the snow intercept. Table A1 lists the specific values 

767 used for the snow mapping. 
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Table 1. Microphysical processes modified or added to the original (i.e., Tao and Simpson 1993; 
Tao et al. 2003) Goddard IM Rutledge and Hobbs-based 3ICE-graupel bulk microphysics 
scheme (updated from L2011). Current changes associated with the new single-moment 4ICE 
scheme are shown in italics. New hail processes are in bold italics. “f( )” indicates “function 
of’. Esi, Egc, and Esc are the collection efficiencies of cloud ice by snow, cloud by graupel and 
cloud by snow, respectively. QcO is the cloud water threshold for snow riming, QiO the cloud ice 
threshold for snow autoconversion, ssi the supersaturation percentage with respect to ice, 
RH/RHice the relative humidity for water/ice, Dwv the water vapor diffusivity, Vs/g the 
snow/graupel fall velocity, Bh,i the immersion mode ice nucleating efficiency, IN the ice number 
concentration, and Tair the air temperature. Qr, Qi, Qs, Qg, and Qh are the rain, cloud ice, snow, 
graupel, and hail mixing ratios, respectively. The process nomenclature essentially follows Ein 
et al. (1983) and Rutledge and Hobbs (1983, 1984). ** =0 if Dhacw + Dhacr < 0.95 * Phwet. 


Process 


Original Modifications 


Referenc e(s)/N otes 


Psaut 


Efficiency Efficiency=0.15, QiO changed from 
f(Tair) g/g to g/m'^, time scale reduced from 

1000 to 300 s, QiO lowered to 
0.4 gm~^, efjiciency=0.25 


Psaci 


Esi = 0.1 Esi f(snow diameter), maximum See snow size mapping in 
Esi = 0.25, maximum Esi=0. 70, Eig. 1 
=0 when Qs=0 
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Praci/Piacr 


Accounts for addition of cloud 


Cloud ice fall speed follows 
Hong et al. (2004) 


Pracw 


Psfi 


Psfw 


Dgacs/Dgaci 


Dgacw 


Dgacr/Pracg 


ice fall speed, =0 when Qr=0, 
becomes hail not graupel 

=0 when Qr=0 

Independent Depends on RH, accounts for Meyers et al. (1992); 


ofRH 

cloud ice size via Meyers IN, 
whieh is a f(ssi), added Dwv 
correction factor & Qi threshold 

Krueger et al. (1995) 


Added Dwv correction factor 



Turned off 

See Lang et al. (2007) 

Ege=1.0 

Ego is f(graupel diameter), 
maximum Ego = 0.65, 

=0 when Qg=0 

See graupel size mapping in 
Eig. 1; Khain et al. (2001) 


=0 when Qg=0 or Qr=0 
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Psacw/Pwacs 


Pracs/Psacr 

Rime 

Splintering 


Pidw/Pidep 


Pint 


Immersion 

Freezing 


Esc=1.0, Esc=0.45, Qe0=1.0 g/kg, 

Qe0=0.5 g/kg Qc0=0.5 g/kg, =0 if Qs=0 


=0 if Qr or Qs=0, =hail not graupel 

None Added and applied to 

Psacw/Pgaew, not f(Vs/g) 
or f(cloud size), added for 

Dhacw 


Based on Based on Meyers IN, which is 
Eletcher a f(ssi), added Dwv 

correction factor to Pidep 

Based on Based on Meyers IN, which is 
Eletcher a f(ssi), previous ice 

concentration checked 

None Added based on Diehl 


Lang et al. (2007); 
Morrison and Grabowski 
(2008) 


Hallet and Mossop (1974); 
f(Tair) and splinter mass 
follow Perrier (1994) 

Eletcher (1962); 

Meyers et al. (1992) 

Eletcher (1962); 

Meyers et al. (1992) 

Diehl and Wurzler (2004); 
Diehl et al. (2006), assumes 
Bh,i= 1.01 e-2 for pollen 
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Contact 


Contaet None 

Added based on Cotton Cotton et al. (1986); 

Nueleation 

and Pruppaeher for Pruppaeher and Klett (1980), 

Brownian diffusion only 500 aetive nuelei per ee 

with radii of 0. 1 mierons 

Saturation Sequential 

Modified sequential, iterative, Tao et al. (2003) 


Adjustment based on Tao ssi up to 10%, ssi up to 15% Reisner (personal) 



(T<-4rC) and 20%, (T>-38°C), 
no evaporation ifW> -0.1 m/s, 
no sublimation if RHice > 70% 
or W > 0 m/s 

Psdep/Pgdep 

=0 if Qs/Qg=0 

Snow/Graupel None 

Allowed if outside eloud and air 

Sublimation 

subsaturated, allowed if air subsaturated 

Pvapg/Pvaph 

convert Qg/Qh to Qs via deposition 
f Qc < l.e- 5 g/g,f(Tair) 

Snow/ Based on 

Based on intereepts mapped aeeording 

Graupel fixed 

to snow/graupel mass and Tair, revised, 
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size 


intercepts 


greater aggregation effect for snow 


Snow density 

Graupel 

density 

Cloud ice 
fall Speed 

Ern 


Pmlts/Pmltg 

Whaci** 

Whacs** 

Whacg** 

Dhacw 


=0.1 g cm'^ =0.05 g ,f (snow size) 

=0.4 g cm'^ =0.3 g crn^ if Qg < 2.0 g m~^ 

=0.5 g cm~^ if Qg > 2.0 g 

None or Based on Hong, included in 
Starr and Cox all sweep volumes 

Added rain evaporation correction 
via tnw, max correction=1.30, 
Added Dwv correction factor 

Added Dwv correction factor 

Added from Lin, =0 if Qh=0, 

Added from Lin, =0 if Qh/Qs=0 

Follows Whacs, =0 if Qh/Qg=0 

Added from Lin, =0 if Qh=0 


Brandes et al. (2007) 

Brown and Swann (1997); 
Straka and Mansell 2005 

Hong et al. (2004); 

Starr and Cox (1985) 

Li et al. (2009) 


Lin et al. (1983) 
Lin et al. (1983) 


Lin et al. (1983) 
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Primh 

convert Qh to Qg via riming, 
f(Phwet, Tair) 

see text on Milbrandt 
and Morrison (2013) 

Dhacr 

Added from Lin, =0 if Qh/Qr=0 

Lin et al. (1983) 

Phwet 

Added from Lin, =0 if Qh=0 

Lin et al. (1983) 

Phfr 

Follows Pgfr but frozen rain 
= hail not graupel, =0 if Qr=0 

Rutledge and Hobbs (1984) 

Phdep 

Added from Lin, =0 if Qh=0 

Lin et al. (1983) 


Phmlt/ Whacr Added from Lin 


Lin et al. (1983) 
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Table 2. Numerical experiments performed for both the 20 May 2012 MC3E and 23 February 
1999 TRMM LB A cases using various options of the Goddard IM bulk microphysics scheme 
where Noh is the hail intercept. 


Experiment Description 


Reference(s)/Notes 


3ice0 3ICE graupel (original) 


Rutledge and Hobbs (1983, 1984) 


3icel Original + no graupel dry collection Lang et al. (2007) 


3ice3 Original + no graupel dry collection 

+ snow/graupel size mapping 


Lang et al. (2011) 


4ice sml New 4ICE with smaller hail 


Noh = 0.0200 cm'' 


4iceb sml New 4ICE with smaller hail 


+ rain evaporation correction 


4iceb med New 4ICE with medium hail 


+ rain evaporation correction 


Noh = 0.0200 cm'' 


Noh = 0.0020 cm' 


-4 


4iceb hg New 4ICE with larger hail 

+ rain evaporation correction 


Noh = 0.0002 cm'' 
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Table Al. Values of snow intercept mapping parameters used to obtain the characteristic snow 
size mapping shown in Fig le. 


Snow Mapping Parameter Cold Value Warm Value Value 


Tcold 



Twarm 



sno\ 

1.0 g m'^ 

0.0 g m 

dsnoX 

4.0 g m'^ 

1.0 g m 

S' exp 1 

1.1 

0.6 


-3 

-3 


5lim 


V, 


sml 


''base 


St exp 1 
tns 

Tslopes (cm"^ C'*) 


-25 “C 
O^C 


0.8 

0.97 


0.040 g m'^ 
0.5 

0.10 cm'"^ 
0.1842 
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FIGURE CAPTIONS 


Figure I. Characteristic sizes (inverse of the slope parameter) of precipitation ice particle 
distributions (inverse exponential) as a function of precipitation ice content and temperature for 
(a) snow in the original Rutledge and Hobbs (1983,1984)-based Goddard scheme, (b) graupel in 
the original Goddard scheme, (c) snow in the modified Goddard SICE scheme, (d) graupel in the 
modified Goddard SICE scheme, (e) snow in the new Goddard 4ICE scheme, (f) graupel in the 
new Goddard 4ICE scheme, and (g) snow density in the new Goddard 4ICE scheme, (a) - (d) 
adapted from E20 1 1 . 

Figure 2. Horizontal cross sections of radar reflectivity for the 20 May 2011 MCSE case (a) 
observed by the NEXRAD Doppler radar network at 10:S0 UTC over north-central Oklahoma 
(figure obtained from the National Mosaic and Next Generation Quantitative Precipitation 
Estimation) and (b) simulated using the new 4ICE scheme with smaller hail and bin rain 
evaporation correction (4iceb sml) at a simulation time of 85.5 h (13:30 UTC). The vertical east- 
west cross section of radar reflectivity shown in (c) was taken through the center of the domain 
from the same simulation and time as (b). 

Figure 3. Time-height cross sections of maximum radar reflectivity for the 20 May 2011 MC3E 
case (a) observed by NEXRAD Doppler radar and simulated using the (b) original 3ICE, (c) 
level 1 improved 3ICE, (d) level 3 improved 3ICE, (e) new 4ICE with smaller hail, (f) new 4ICE 
with smaller hail and bin rain evaporation correction, (g) new 4ICE with moderate hail and bin 
rain evaporation correction, and (h) new 4ICE with larger hail and bin rain evaporation 
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correction Goddard microphysics scheme. Right axes are heights in km, while horizontal dashed 
lines show the level of indicated environmental temperatures in degrees C. Times range from 00 
UTC 20 May to 00 UTC 2 1 May 2011. 

Figure 4. Same as Fig. 2c except showing vertical cross sections of simulated (a) vertical 
velocities (b) graupel fall speeds, and (c) hail fall speeds. 

Figure 5. Vertical profiles of maximum radar refieetivity for the 20 May 2011 MC3E case 
extracted between 06 and 12 UTC Ifom Doppler radar observations and between 09 and 15 UTC 
Ifomthe three Goddard SICE simulations and four Goddard 4ICE simulations shown in Figure 2. 

Figure 6. Radar refieetivity CFADs for the 20 May 2011 MC3E case constructed Ifom (a) 
NEXRAD Doppler radar observations and simulations using the (b) original SICE, (c) level 1 
improved SICE, (d) level 3 improved SICE, (e) new 4ICE with smaller hail, (f) new 4ICE with 
smaller hail and bin rain evaporation eorrection, (g) new 4ICE with moderate hail and bin rain 
evaporation eorrection, and (h) new 4ICE with larger hail and bin rain evaporation eorrection 
Goddard microphysics scheme. Heavy thick lines in (b) - (h) show the edges of the core 
observed Ifequeney probabilities [i.e., the 5 % eontours shown in (a)] and the outer limits of the 
observed Ifequeney distributions [i.e., the 0 % eontours shown in (a)]. Right axes are heights in 
km, while horizontal dashed lines show the level of indicated environmental temperatures in 
degrees C. 
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Figure 7. Vertical profiles of PDF matching scores (i.e., the amount of overlap between the 
simulated and observed PDF at each level) for the 20 May 2011 MC3E simulations using the (b) 
original SICE, (c) level 1 improved SICE, (d) level S improved SICE, (e) new 4ICE with smaller 
hail, (f) new 4ICE with smaller hail and bin rain evaporation correction, (g) new 4ICE with 
moderate hail and bin rain evaporation correction, and (h) new 4ICE with larger hail and bin rain 
evaporation correction Goddard microphysics scheme. 

Eigure 8. Horizontal cross sections of radar reflectivity for the 23 Eebruary 1999 LBA case (a) 
observed by the S-Pol radar at 20:50 UTC over Amazonia overlaid with storm-relative winds 
Ifom a dual-Doppler wind analysis and the track of the University of North Dakota Citation 
aircraft (figure adapted Ifom http://radarmet.atmos.colostate.edu/lba_trmm/23feblba_cappi.html 
and Lang et al. 2007) and (b) simulated using the new 4ICE scheme with smaller hail and bin 
rain evaporation correction (4iceb sml) at a simulation time of 330 minutes (21:00 UTC). The 
vertical east-west cross section of radar reflectivity shown in (c) was taken through the center of 
the domain Ifom the same simulation and time as (b). The solid rectangle and dashed box shown 
in (b) denote the north-south oriented rectangular patch of higher sensible/lower latent heat 
fluxes (Ji Parana) imposed to initiate convection and the analysis domain, respectively. 

Eigure 9. Time-height cross sections of maximum radar reflectivity for the 23 Eebruary 1999 
LBA case (a) observed by the S-pol ground-based radar and simulated using the (b) original 
3ICE, (c) level 1 improved 3ICE, (d) level 3 improved 3ICE, (e) new 4ICE with smaller hail, (f) 
new 4ICE with smaller hail and bin rain evaporation correction, (g) new 4ICE with moderate hail 
and bin rain evaporation correction, and (h) new 4ICE with larger hail and bin rain evaporation 
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correction Goddard microphysics scheme. Right axes are heights in km, while horizontal dashed 
lines show the level of indicated environmental temperatures in degrees C. Model data were 
taken from a 64 km x 64 km subdomain. Blaek and gray labels at the bottom of (a) are the UTC 
and approximate matching times, respectively. 

Figure 10. Same as Fig. 8c except showing vertieal cross sections of simulated (a) vertical 
velocities (b) graupel fall speeds, and (c) hail fall speeds. 

Figure 11. Vertical profiles of the maximum radar reflectivity for the 23 February 1999 LBA 
ease extraeted from the S-pol radar observations and the last 60 minutes of the three Goddard 
SICE simulations and four Goddard 4ICE simulations shown in Eigure 6. Model data were taken 
from a 64 km x 64 km subdomain. 

Eigure 12. Radar reflectivity CFADs for the 23 Eebruary 1999 LBA case constructed from (a) S- 
pol radar observations and the final 60 minutes of the simulations using the (b) original SICE, (c) 
level 1 improved SICE, (d) level 3 improved SICE, (e) new 4ICE with smaller hail, (f) new 4ICE 
with smaller hail and bin rain evaporation correction, (g) new 4ICE with moderate hail and bin 
rain evaporation correction, and (h) new 4ICE with larger hail and bin rain evaporation 
correction Goddard microphysics scheme. Heavy thick lines in (b) - (h) show the edges of the 
core observed frequency probabilities [i.e., the 5 % contours shown in (a)] and the outer limits of 
the observed frequency distributions [i.e., the 0 % contours shown in (a)]. Right axes are heights 
in km, while horizontal dashed lines show the level of indicated environmental temperatures in 
degrees C. Model data were taken from a 64 km x 64 km subdomain. 


56 


1480 

1481 

1482 

1483 

1484 

1485 

1486 

1487 

1488 

1489 

1490 

1491 

1492 

1493 

1494 

1495 

1496 

1497 

1498 

1499 

1500 

1501 

1502 


Figure 13. Vertical profiles of PDF matching scores for the 23 February 1999 LBA simulations 
over the final 60 minutes using the (b) original 31CE, (c) level 1 improved 31CE, (d) level 3 
improved 31CE, (e) new 41CE with smaller hail, (f) new 41CE with smaller hail and bin rain 
evaporation correction, (g) new 41CE with moderate hail and bin rain evaporation correction, and 
(h) new 41CE with larger hail and bin rain evaporation correction Goddard microphysics scheme. 
Model data were taken Ifom a 64 km x 64 km subdomain. 

Eigure 14. Instantaneous surface rainfall rates corresponding to the horizontal radar reflectivity 
cross sections shown for the 20 May MC3E case in Fig. 2b (a) and the 23 February EBA case in 
Eig. 8b (b). (c) surface rainfall histograms observed by the Doppler radar network around the 
MC3E sounding array Ifom 06-12 UTC and simulated with Goddard microphysics Ifom 09-15 
UTC for the 20 May MC3E case, (d) surface rainfall histograms derived from ground-based 
radar observations collected Ifom 2002-2130 UTC and simulated over the final 60 minutes of 
simulation time over a 64 km x 64 km subdomain (shown by the dashed square in panel b) for 
the 23 Eebruary 1999 case using the Goddard microphysics schemes. 

Eigure 15. Distribution of surface cold pool intensities for the 20 May MC3E case for the 
smaller hail runs with and without the bin rain evaporation correction, 4iceb sml and 4ice sml, 
respectively. Intensities are shown in terms of the surface potential temperature deviations (K) 
over the 09 to 15 UTC analysis period for regions where the lowest level rain mixing ratio 
exceeds 0.1 g m'^. 
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1546 Figure 1. Characteristic sizes (inverse of the slope parameter) of precipitation ice particle 

1 547 distributions (inverse exponential) as a function of precipitation ice content and temperature for 

1548 (a) snow in the original Rutledge and Hobbs (1983,1984)-based Goddard scheme, (b) graupel in 

1 549 the original Goddard scheme, (c) snow in the modified Goddard SICE scheme, (d) graupel in the 

1 550 modified Goddard SICE scheme, (e) snow in the new Goddard 4ICE scheme, (f) graupel in the 

1551 new Goddard 4ICE scheme, and (g) snow density in the new Goddard 4ICE scheme, (a) - (d) 

1 552 adapted from E20 1 1 . 
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Figure 2. Horizontal cross sections of radar reflectivity for the 20 May 2011 MC3E case (a) 
observed by the NEXRAD Doppler radar network at 10:30 UTC over north-central Oklahoma 
(figure obtained Ifom the National Mosaic and Next Generation Quantitative Precipitation 
Estimation) and (b) simulated using the new 4ICE scheme with smaller hail and bin rain 
evaporation correction (4iceb sml) at a simulation time of 85.5 h (13:30 UTC). The vertical east- 
west cross section of radar reflectivity shown in (c) was taken through the center of the domain 
Ifom the same simulation and time as (b). 
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Figure 3. Time-height cross sections of maximum radar reflectivity for the 20 May 2011 MC3E 
case (a) observed by NEXRAD Doppler radar and simulated using the (b) original 31CE, (c) 
level 1 improved 31CE, (d) level 3 improved 31CE, (e) new 41CE with smaller hail, (f) new 41CE 
with smaller hail and bin rain evaporation correction, (g) new 4ICE with moderate hail and bin 
rain evaporation correction, and (h) new 41CE with larger hail and bin rain evaporation 
correction Goddard microphysics scheme. Right axes are heights in km, while horizontal dashed 
lines show the level of indicated environmental temperatures in degrees C. Times range from 00 
UTC 20 May to 00 UTC 2 1 May 2011. 
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Figure 4. Same as Fig. 2e except showing vertical cross sections of simulated (a) vertical 
velocities (b) graupel fall speeds, and (c) hail fall speeds. 
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Figure 5. Vertical profiles of maximum radar reflectivity for the 20 May 2011 MC3E case 
extracted between 06 and 12 UTC from Doppler radar observations and between 09 and 15 UTC 
from the three Goddard SICE simulations and four Goddard 4ICE simulations shown in Figure 2. 
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Figure 6. Radar reflectivity CFADs for the 20 May 2011 MC3E case constructed from (a) 
NEXRAD Doppler radar observations and simulations using the (b) original 31CE, (c) level 1 
improved 31CE, (d) level 3 improved 31CE, (e) new 41CE with smaller hail, (f) new 41CE with 
smaller hail and bin rain evaporation correction, (g) new 41CE with moderate hail and bin rain 
evaporation correction, and (h) new 41CE with larger hail and bin rain evaporation correction 
Goddard microphysics scheme. Heavy thick lines in (b) - (h) show the edges of the core 
observed frequency probabilities [i.e., the 5 % contours shown in (a)] and the outer limits of the 
observed frequency distributions [i.e., the 0 % contours shown in (a)]. Right axes are heights in 
km, while horizontal dashed lines show the level of indicated environmental temperatures in 
degrees C. 
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1625 Figure 7. Vertical profiles of PDF matching scores (i.e., the amount of overlap between the 

1 626 simulated and observed PDF at each level) for the 20 May 2011 MC3E simulations using the (b) 

1627 original SICE, (c) level 1 improved SICE, (d) level 3 improved SICE, (e) new 4ICE with smaller 

1628 hail, (f) new 4ICE with smaller hail and bin rain evaporation correction, (g) new 4ICE with 

1 629 moderate hail and bin rain evaporation correction, and (h) new 4ICE with larger hail and bin rain 

1 630 evaporation correction Goddard microphysics scheme. 
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Figure 8. Horizontal cross sections of radar reflectivity for the 23 February 1999 LBA case (a) 
observed by the S-Pol radar at 20:50 UTC over Amazonia overlaid with storm-relative winds 
from a dual-Doppler wind analysis and the track of the University of North Dakota Citation 
aircraft (figure adapted from http://radarmet.atmos.colostate.edu/lba_trmm/23feblba_cappi.html 
and Lang et al. 2007) and (b) simulated using the new 4ICE scheme with smaller hail and bin 
rain evaporation correction (4iceb sml) at a simulation time of 330 minutes (21:00 UTC). The 
vertical east -west cross section of radar reflectivity shown in (c) was taken through the center of 
the domain from the same simulation and time as (b). The solid rectangle and dashed box shown 
in (b) denote the north-south oriented rectangular patch of higher sensible/lower latent heat 
fluxes (Ji Parana) imposed to initiate convection and the analysis domain, respectively. 
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Figure 9. Time-height cross sections of maximum radar reflectivity for the 23 February 1999 
LBA case (a) observed by the S-pol ground-based radar and simulated using the (b) original 
31CE, (c) level 1 improved 31CE, (d) level 3 improved 31CE, (e) new 41CE with smaller hail, (f) 
new 41CE with smaller hail and bin rain evaporation correction, (g) new 41CE with moderate hail 
and bin rain evaporation correction, and (h) new 41CE with larger hail and bin rain evaporation 
correction Goddard microphysics scheme. Right axes are heights in km, while horizontal dashed 
lines show the level of indicated environmental temperatures in degrees C. Model data were 
taken from a 64 km x 64 km subdomain. Black and gray labels at the bottom of (a) are the UTC 
and approximate matching times, respectively. 
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Figure 10. Same as Fig. 8e except showing vertical cross sections of simulated (a) vertical 
velocities (b) graupel fall speeds, and (c) hail fall speeds. 
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1681 Figure 11. Vertical profiles of the maximum radar reflectivity for the 23 February 1999 LBA 

1682 case extracted Ifom the S-pol radar observations and the last 60 minutes of the three Goddard 

1 683 31CE simulations and four Goddard 41CE simulations shown in Eigure 6. Model data were taken 

1684 Ifom a 64 km x 64 km subdomain. 
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Figure 12. Radar reflectivity CFADs for the 23 February 1999 LBA case constructed from (a) S- 
pol radar observations and the final 60 minutes of the simulations using the (b) original 31CE, (c) 
level 1 improved 31CE, (d) level 3 improved 31CE, (e) new 41CE with smaller hail, (f) new 41CE 
with smaller hail and bin rain evaporation correction, (g) new 41CE with moderate hail and bin 
rain evaporation correction, and (h) new 41CE with larger hail and bin rain evaporation 
correction Goddard microphysics scheme. Heavy thick lines in (b) - (h) show the edges of the 
core observed frequency probabilities [i.e., the 5 % contours shown in (a)] and the outer limits of 
the observed frequency distributions [i.e., the 0 % contours shown in (a)]. Right axes are heights 
in km, while horizontal dashed lines show the level of indicated environmental temperatures in 
degrees C. Model data were taken from a 64 km x 64 km subdomain. 
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1747 Figure 13. Vertical profiles of PDF matching scores for the 23 February 1999 LBA simulations 

1748 over the final 60 minutes using the (b) original 31CE, (c) level 1 improved 31CE, (d) level 3 

1749 improved 31CE, (e) new 41CE with smaller hail, (f) new 41CE with smaller hail and bin rain 

1750 evaporation correction, (g) new 41CE with moderate hail and bin rain evaporation correction, and 

1751 (h) new 41CE with larger hail and bin rain evaporation correction Goddard microphysics scheme. 

1752 Model data were taken Ifom a 64 km x 64 km subdomain. 
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1756 Figure 14. Instantaneous surface rainfall rates corresponding to the horizontal radar reflectivity 

1757 cross sections shown for the 20 May MC3E case in Fig. 2b (a) and the 23 February LBA case in 

1758 Fig. 8b (b). (c) surface rainfall histograms observed by the Doppler radar network around the 

1759 MC3E sounding array from 06-12 UTC and simulated with Goddard microphysics from 09-15 

1760 UTC for the 20 May MC3E case, (d) surface rainfall histograms derived from ground-based 

1761 radar observations collected from 2002-2130 UTC and simulated over the final 60 minutes of 

1762 simulation time over a 64 km x 64 km subdomain (shown by the dashed square in panel b) for 

1763 the 23 Eebruary 1999 case using the Goddard microphysics schemes. 
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1767 Figure 15. Distribution of surface cold pool intensities for the 20 May MC3E case for the 

1768 smaller hail runs with and without the bin rain evaporation correction, 4iceb sml and dice sml, 

1769 respectively. Intensities are shown in terms of the surface potential temperature deviations (K) 

1770 over the 09 to 15 UTC analysis period for regions where the lowest level rain mixing ratio 

1771 exceeds 0. 1 g m'^. 
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